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Abstract. Forming stars emit a significant amount of radiation into their natal environment. 
While the importance of radiation feedback from high-mass stars is widely accepted, radiation 
has generally been ignored in simulations of low-mass star formation. I use ORION, an adaptive 
mesh refinement (AMR) three-dimensional gravito-radiation-hydrodynamics code, to model low- 
mass star formation in a turbulent molecular cloud. I demonstrate that including radiation 
feedback has a profound effect on fragmentation and protostellar multiplicity. Although heating 
is mainly confined within the core envelope, it is sufficient to suppress disk fragmentation that 
would otherwise result in low-mass companions or brown dwarfs. As a consequence, turbulent 
fragmentation, not disk fragmentation, is likely the origin of low-mass binaries. 
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1. Introduction 

The origin of stellar multiplicity remains an unsolved problem in star formation. The 
dense conditions in molecular cloud cores and dim luminosities of young protostars make 
estimations of the initial multiplicity distribution challenging (Duchene et al. 2007). How- 
ever, the present-day multiplicity can be observed among field stars, where the likelihood 
of companions is strongly correlated with the primary stellar mass. Nearly all O and B 
stars are found in binaries or multiple systems, while only ^20% of M stars have compan- 
ions (Lada 2006). Successful simulations and theories of star formation must be able to 
predict the multiplicity fraction and explain why it depends so strongly on stellar mass. 

A number of mechanisms have been proposed of which two appear to have the most 
potential for producing the observed number of multiple systems (see Tohline 2002 for a 
review). First, gravitational instability within a protostellar accretion disk may produce 
companions within a few 100 AU (Adams et al. 1989; Bonnell et al. 1994). Repeated 
fragmentation over the disk lifetime may generate numerous companions. Second, per- 
turbations within a turbulent core may seed additional fragmentation on scales of ~ 
0.001-0.1 pc (Fisher 2004; Goodwin et al. 2004). This must occur within the first ~ 0.5 
Myr of collapse, when the core still contains at least 0.1M Q . 

In this paper, we discuss the effect of radiation feedback on early fragmentation and 
stellar multiplicity using 3D adaptive mesh refinement (AMR) simulations of turbulent 
molecular clouds. We compare the cases with and without radiation from forming stars. 
We describe the simulations in §2, present results in §3, and conclude in §4. 



2. Simulations 

The ORION code solves the equations of compressible gas dynamics, Poisson equation, 
and radiation energy equation in the flux-limited diffusion approximation (Krumholz et 



al. 2007): 



dt 
d{pv) 
dt 
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^+V.(pv)^0, (2.1) 



V • (pvv) = -VP - pV0, (2.2) 



+ V • [(pe + P)v] = pvV0 - KR p(4irB - cE), (2.3) 

— - V • ( VP) = k pP {AttB - cE) + V L„<S(x - x n ), (2.4) 

dt re R p ^ 

V 2 = 4^G[p + ^m„<5(x-x n )], (2.5) 

where p, P, v and e are the fluid density, pressure, velocity, and specific kinetic plus 
internal energy of the gas, <f> is the gravitational potential, m n , x n , and L n are the mass, 
position, and luminosity of the n th star, E is the radiation energy density, and kr and 
rep are the Rosseland and Planck dust opacities. For comparison, a second calculation 
closes the equations with a barotropic equation of state (EOS) in lieu of equation 12.41 



P = pc 2 s + {£) P^l (2-6) 

where c s — {k^T / /Lt) 1 / 2 is the sound speed, 7 = 5/3, the average molecular weight 
jtx = 2.33mH, and the critical density, p c = 2x 10~ 13 g cm~ 3 . 

The calculations insert Lagrangian sink particles in regions exceeding the Jeans density 
on the maximum AMR level (Krumholz et al. 2004). In the radiation simulation, particles 
have a sub-grid model based upon McKee & Tan (2003) that includes the accretion 
energy, Kelvin- Helmholtz contraction, and nuclear burning (Offher ct al. 2009). 

The calculations have a Mach number of A^3D=6.6, domain size L = 0.65 pc, and mass 
M=185 Mq, which correspond to an approximately virialized cloud. We adopt periodic 
boundary conditions for the gas and Marshak boundary conditions for the radiation 
field, which allows the cloud to cool. We use a 256 3 base grid with 4 levels of grid 
refinement, where AX4 — 32 AU. High- resolution convergence tests are discussed in 
detail in Offher et al. (2009). We drive the boxes for three crossing times using random 
velocity perturbations with wavenumbers 1 ^ k ^ 2 after which self-gravity is turned 
on. The initial gas temperature is 10 K. Since the gas cools efficiently during the driving 
phase, the radiation calculation remains nearly isothermal. 



3. Results 

3.1. Fragmentation 

At the end of a freefall time, the two calculations have very different temperature distri- 
butions. Temperatures in the EOS calculation do not exceed 15 K, while temperatures 
in the radiation calculation reach ~ 100 K. In the latter, because the heated cores are 
turbulent, a range of densities is heated to various temperatures such that gas temper- 
ature is not a single valued function of gas density. The heating also varies with the 
number of stars and their instantaneous accretion rate, such that any EOS fit to a mean 
temperature-density curve would be a generally poor approximation over the course of 
the simulation. 

As a consequence of the feedback, the calculation with radiation has a fundamentally 
different mass distribution (Figure [T]A_). The protostellar heating is sufficient to raise the 
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Figure 1. Left: Distribution of star masses produced from simulations with (dashed) and with- 
out (solid) radiative feedback. Right: Pair separation as a function of time in 1 kyr bins for all 
particle pairs from Offner et al. (2010). The dashed line at 500 AU indicates a rough boundary 
between the disk and core scales. The large majority of pairs have separations > 0.1 pc and are 
not shown. The large (red) symbols indicate the first time bin. 

temperature of the gas significantly within a few hundred AU and suppress fragmentation 
in accretion disks that would otherwise be unstable. As a result, the radiative calcula- 
tion has far fewer brown dwarfs and, despite small number statistics, its mass distribu- 
tion more closely resembles the stellar initial mass function. Bate (2009) has previously 
demonstrated similar results using smoothed-particle hydrodynamics simulations. 

Although the heating is efficient at small scales, it is limited to the parent core and does 
not inhibit fragmentation of other cores or even wide fragmentation within the same core. 
Figure QJ3 shows the separations of all star pairs as a function of time in the calculation 
with radiation feedback. The plot is restricted to separations of 0.1 pc or less to highlight 
stars forming within the same core that may comprise a binary. The plot shows that 
fragmentation on scales of ~ 2000 AU is not suppressed by heating. The one instance 
of fragmentation with d < 500 AU is actually filament rather than disk fragmentation. 
This suggests that turbulent core fragmentation remains a viable mechanism for binary 
formation of low-mass stars, while disk fragmentation is much less likely. 

3.2. Disk Analysis 

Although the disks in the simulations are not well resolved, it is still possible to use their 
mean properties to draw robust conclusions about the protostellar multiplicity. Kratter 
et al. (2010) define a two-dimensional parameter space for characterizing accretion and 
disk stability: 



M {n G 



r = 



(3-1) 



The thermal parameter, £, compares the core sound speed to the disk sound speed, c Sj <j, 
where M- m is the infall mass accretion rate. A collapsing isothermal sphere has £ ~ 1 (Shu 
1977). For £ > 1, a disk will be unable to efficiently process accreting material and will 
eventually fragment. The rotational parameter, T, compares the disk orbital time to the 
gas infall time, where M*d is the total mass in the star-disk system, f2k,in is the Keplerian 
angular velocity at the circularization radius of the infall, and (j) in is the specific angular 
momentum. For large r (r ~ 0.1) the disk mass changes quickly over an orbital time. 
In Figure[2]we estimate these parameters for the disks in each calculation. Our analysis 
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Figure 2. The values of V and £ at 1 freefall time for the protostellar disks in each simulation, 
where the Non-Radiative Transfer (NRT) cases are denoted by boxes. The diagonal line indicates 
the boundary between stable and unstable disks found by Kratter et al. (2010). 

confirms that disks in the radiative calculation are stable, where the binaries to the left 
of the line are in fact products of turbulent rather than disk fragmentation. Without 
radiation feedback large unstable disks yield high multiplicity systems, which fall to the 
right of the line as expected. Although some single systems exist within the unstable 
regime these tend to have low Q values and often have previously fragmented. 

4. Conclusions 

Radiation feedback from low-mass stars is important in shaping the stellar mass dis- 
tribution. Heating works to stabilize protostellar disks and suppress fragmentation that 
would otherwise over-produce brown dwarfs. However, fragmentation of the parent core 
may still occur on thousand AU scales, suggesting that turbulent core fragmentation, not 
disk fragmentation, is the most likely origin of low-mass binaries. 
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